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Abstract 


This thesis develops a method to model the acoustic field generated by a 
monopole source placed in a moving rectangular duct. The walls of the duct are assumed 
to be infinitesimally thin and the source is placed at the center of the duct. The total 
acoustic pressure is written in terms of the free-space pressure, or incident pressure, and 
the scattered pressure. The scattered pressure is the augmentation to the incident pressure 
due to the presence of the duct. It satisfies a homogeneous wave equation and is 
discontinuous across the duct walls. Utilizing an integral representation of the scattered 
pressure, a set of singular boundary integral equations governing the unknown jump in 
scattered pressure is derived. This equation is solved by the method of collocation after 
representing the jump in pressure as a double series of shape functions. The solution 
obtained is then substituted back into the integral representation to determine the 
scattered pressure, and the total acoustic pressure at any point in the field. A few 
examples are included to illustrate the influence of various geometric and kinematic 
parameters on the radiated sound field. 
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1. Introduction 


The understanding and accurate prediction of sound radiation associated with 
aircraft applications is of current interest. Recently, a great deal of attention has been 
directed to modeling the ultrahigh by-pass ratio turbo-fan engine (ducted propfan) due to 
its efficiency and shrouded propeller design. This engine has become attractive as an 
efficient component of future commercial transport aircraft. It is known that an 
unshrouded propeller generates an acoustic field which tends to radiate in the lateral 
direction. By housing the propeller in a shroud or duct (as in the ducted propfan), 
potential benefits of noise reduction exist since the duct provides a shield in the primary 
radiation direction. The ability to accurately model this propeller noise in both the ducted 
and unducted cases is of great importance. Modeling techniques such as computational 
fluid dynamics (CFD), finite element methods (FEM) and boundary integral techniques 
have been utilized to predict the acoustic benefits of the ducted propfan. Eversman [1] 
developed a finite element model for the generation, propagation and radiation of noise of 
a ducted fan. He also constructed a ffee-field propeller model compatible with the finite 
element formulation and conducted noise studies as presented in ref. [1]. Radiated field 
results for both the ducted and unducted propeller were obtained. Lan [2] and Buhler [3] 
developed separate prediction methods based on a boundary integral technique for the 
acoustic field generated by a propeller within a circular, rigid duct. 

The presence of a duct or shroud was shown in all above studies to reduce the 
level of noise radiated from a propeller, at least in the primary radiation direction. With 
the presence of the shroud, an acoustically treated duct wall can also be util' zed to furthei 
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increase the level of noise reduction obtained. The work of Kosanchick [4] was an 
extension of the work presented in refs. [2,3] and included an acoustically treated duct 
wall. Dunn, Tweed and Farassat [5] also presented a similar model of the acoustic field 
generated by a ducted propfan. Their boundary integral equation technique modeled 
acoustically treated walls, but only rigid wall results were presented in ref. [5]. The 
models in refs. [2-5] required fewer overall computations than the finite element 
technique since the boundary integral method only requires dealing with a surface 
integration along the duct as compared to finite element computations over the entire 
region of interest. However, the boundary integral technique cannot be applied if 
complex geometrical models or complex mean flows are involved. 

As it is of great importance to develop accurate prediction methods of sound 
radiation from a ducted propeller, it is also of importance to understand the relationship 
between the duct geometry and source type in regard to noise reduction. The analysis in 
refs. [2-4] was modified and applied to a point source rather them a propeller model by 
Myers [6]. With a simpler source model, comparisons could be made between the 

boundary integral technique and the CFD approach. Ozyoriik and Long [7] utilized a 
finite-difference approach in solving this particular problem, and their results were found 
to be in agreement with those of ref. [6], 

The studies previously discussed in refs. [1-6] dealt with only an axisymmetric 
duct configuration. However, rectangular duct configurations are often utilized for 
experimental inlet research. Thus, it is the purpose of the current work to develop a 
boundary integral technique similar to those in refs. [2-6] to describe radiation from a 
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monopole or point source placed within a rectangular duct. This particular duct geometry 
lacks symmetry and requires significantly greater computational effort than the problems 
previously discussed in refs. [2-4,6]. The source was restricted to be a monopole so that 
at least the incident field would be symmetric in form. The duct walls are assumed to be 
rigid and infinitesimally thin. The current work utilizes a scattering formulation which 
decouples the effects of the duct and source [2-6]. The scattered pressure is 
discontinuous across the walls of the duct. It is governed by a generalized wave equation 
with a source term that is proportional to its unknown jump across the duct walls. An 
integral representation for the solution of this in terms of this jump in scattered pressure 
across the duct surface is obtained. The scattered pressure jumps are then represented by 
two carefully chosen shape function expansions. The free-space pressure is used as the 
input data and the method of collocation is utilized to solve a system of algebraic 
equations for the coefficients of the shape function expansions. From these coefficients, 
the jump in scattered pressure along the duct surface is obtained and substituted into the 
original integral representation to obtain the scattered pressure. The total acoustic 
pressure is obtained through the sum of the free-space or incident pressure and the 

scattered pressure for a particular field point. 

The remainder of this thesis discusses the formulation and validation of the 
boundary integral technique in its application to the rectangular duct with a monopole 
source at its center. In particular, Chapter 2 details the development of the governing 
boundary integral equation. Chapter 3 presents the numerical details necessary to 
accurately obtain the solution to the set of algebraic equations. The impact of the choice 
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of parameters on the solution is shown through a few example problems as discussed in 
Chapter 4. 
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2. Theory 

This chapter presents the development of the governing equations required for the 
prediction of the acoustic field radiated from a source in a moving rectangular duct. 

2.1 Integral Representation 

An acoustic source in an infinitesimally thin, rigid, rectangular duct is examined. 
The rectangle is of length L, width 2b and height 2a as shown in Figure 1 . A body-fixed 
cartesian coordinate system X , placed at the center of the moving duct, and an inertially- 
fixed cartesian coordinate system x , are used to describe the duct. The duct is assumed 
to be moving subsonically at velocity V, in the negative x 3 direction of the inertially-fixed 
frame. The objective is to obtain the acoustic field radiated to free space through the 
open ends of the duct. 

By linearizing the ideal fluid equations of motion it can be shown that acoustic 
wave propagation in isentropic flow with no body forces is governed by the wave 
equation which can be written in the form 

□ 2 p,--T§ L - v 2 p.=q( s - t ) (21) 

c oi 

where p, is the total acoustic pressure, q(x,t)is the noise source, and c is the speed of 
sound in the fluid medium. 

The total acoustic pressure is written in terms of incident and scattered acoustic 
pressure components as 
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where p ( is the free-space pressure due to the monopole source alone and p s is the 
augmentation of the incident pressure due to the presence of the rectangular duct. 
Through the use of this scattering formulation, the effects of the duct and noise source are 
decoupled and can be handled separately. The incident field satisfies 

□ 2 Pi =q(x,t) (2.3) 

It then follows that the scattered field must satisfy the homogeneous wave equation 

□ 2 p a = 0 (2.4) 

The objective is to obtain an integral representation for the scattered field along with the 

appropriate boundary conditions at the duct surface and a radiation condition specifying 

that the scattered field is outgoing in the region exterior to the duct. 

Let f(x)=0 describe the surface of the thin, rigid duct moving in a direction 
tangent to itself where f is defined such that Vf = n and n is the unit outward normal to 
the duct surface. The presence of the solid duct walls gives rise to a discontinuity in 
scattered pressure across the duct surface. Through the use of generalized derivatives and 
eqn.(2.4), it can be shown that the scattered pressure is a solution to the generalized wave 
equation 

□ 2 P S = v[Ap s n8(f)] (2.5) 

where §(•) is the Dirac delta function and Ap s is the jump in scattered pressure across the 
duct surface [8]. The bars over the differential operators in eqn. (2.5) signify generalized 
differentiation. Equation (2.5) is a special form of the well-known Ffowcs Williams- 
Hawkings (FW-H) equation [8,9]. 
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By solving eqn.(2.5) utilizing the free-space Green’s function for the wave 


equation, an integral representation for p s is obtained in the form 


1 3 


Ap s cos0 


r 1 - M, 


dS- 


Ap s cos0 


r 2 1 - M, 


dS 


( 2 . 6 ) 


f = 0 f = 0 

The details of the derivation of eqn. (2.6) are given in number of earlier publications [8] 
and therefore will not be repeated here. In eqn. (2.6), r = |r| = |x - y| is the magnitude of the 
radiation vector representing the distance from a source point at y on f=0 to an observer 


at x ; M r is the component of the surface Mach number in the direction of the radiation 
vector; 0 is the angle between f and the unit normal to the surface of the duct n at the 
source location; dS is the elemental area of the duct surface f^O. The integrands are 
evaluated at the emission time t* which is a solution to the retarded time equation 


There is only one solution t for the retarded time equation since the duct is assumed to 
be moving subsonically through the fluid medium. Physically, r is the emission time of 
signals generated by a source on the duct surface that are received by an observer at 
position x at time t. 

By applying the appropriate boundary conditions to eqn. (2.6), an integral 
equation will be developed from which the unknown jump in scattered pressure Ap s is 
determined over the surface of the duct. Once A p s is known, the scattered acoustic field 
for any position and time (x,t) can be determined using eqn. (2.6) again. 
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The following details of the analysis are simplified if eqn. (2.6) is expressed in 
terms of the body-fixed coordinate system. To relate the observer and source positions in 
the body-fixed reference frame to the observer and source positions in the inertially-fixed 
reference frame, the following transformation is utilized: 

X - S+V,T " ( 2 . 8 ) 

Y-y + Vti, 

Here T s is the unit vector in the 3-direction. 

The analytical development of the boundary integral equation to follow can be 
carried out for an arbitrary incident pressure field p. However, computational 
complexities make it desirable to introduce some symmetry into the problem. Thus, at 
this stage it is assumed that the acoustic source is a monopole and, for further simplicity, 
it is positioned at the origin of the body-fixed coordinate system (i.e., at the center of the 
duct). Among other things, this restriction ensures that the scattered pressure jumps are 
identical on surfaces x 2 = ±b and Xi = ±a . The incident field corresponding to a moving 
monopole source is reviewed in the next section. 

2.2 Incident Monopole Field 

The incident field is the field due to the source alone without the presence of the 
duct. The solution to eqn. (2.3) can be obtained for the monopole source in terms of a 
complex velocity potential that satisfies 

1 ^ (h 

□ 2 <t>j = ~ ~ZT~ V 2 ^ = Ae' ,a>, 8(x 3 + Vt)5(x 2 )5(xi) (2.9) 

c at 
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where A is defined as the source strength. The solution to eqn. (2.9) is well known and 
appears, for example, in Morse and Ingard [10]. In terms of the body-fixed coordinates 
the solution presented in [10] is 


- Ae' i<Dt 

4>(X,t) = — /^ exp| iaj 


4txVb^ 


-mx 3 
V p 


y[K 


( 2 . 10 ) 


where the symbol a, which is used throughout this work, is defined by a = co/cP , and 
p 2 = 1 - m 2 . Use of the relation 

3*1 


Pi = -Po 


dt 


= -p. 

r dt 

a i» ") 
+ v-^- 

* 0 1 
X 

v dt 

X dxj 


( 2 . 11 ) 


yields the incident monopole complex pressure field in the form 


P, : 


-Po cA e~ 
4k p 2 


-MX, 
ice — - — + 


Bo 




MX, 


(Bo)' 


e V p 




( 2 . 12 ) 


where B 0 = — + xf + xl and p 0 is the fluid density. As noted above, the monopole 
P 2 

complex pressure field expressed in eqn. (2.12) will be the only incident field considered 
here. As with all quantities in the following, the physical pressure is taken to be the real 
part of the complex incident pressure. 


2.3 Boundary Integral Formulation 


The analysis involved in formulating the boundary integral equation requires 
lengthy algebraic manipulations of eqn. (2.6) that are simplified if it is written in terms of 
the body-fixed coordinates. For convenience define £, r\ and C, as 
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ri = Y 2 -X 2 , ? = Yi-Xi 


( 2.13) 


r Y 3 ~ X 3 

P 

These abbreviations are utilized frequently throughout this work. The radiation vector, 
which extends from the duct surface f=0 to an observer at x , is 

r = x-y = X- Y - V(t-T)fj (2.14) 

Substitution of eqn. (2.14) into the retarded time equation eqn.(2.7), leads to a quadratic 
equation for the radiation distance at time t\ whose solution is 

= c(t- T *) (2.15) 

The component of the surface Mach number in the direction of the radiation vector is 


(ms + ^W+C 2 ) 


f -Tm^P + M 2 r 

Mr = M/3- — = — 1 

r r 


(2.16) 


where r is the magnitude of the radiation vector. Evaluation of eqn. (2. 1 6) at the emission 
time x* and use of eqn. (2. 1 5) yields 

Pjg+g+g 


I -M r 


(2.17) 


The factor cos0 is given by 


cos0 = ^ ( 2.18) 

r 

Finally, the scattered pressure jump Ap s can be modeled in the form suggested by the 
incident field as 


Ap, ( X. » t) = «( X s ) e'™ 1 (2.19) 

where X s represents surface coordinates. Evaluation of eqn. (2.19) at the emission time 

gives 
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/ \ / \ io>r* 

A Ps(Y,,'t*) = it(Y,)e'“‘eT' e p 


( 2.20) 


Substitution of eqns. (2.17), (2.18) and (2.20) into the two surface integrals representing 
the scattered pressure in eqn. (2.6), then puts the integral representation of p s in the form 


] 8 

47tP,(X,t) = --- 


rc(Y s )n*r e lwt e c e 

r*V^ + r l 2 + ? 2 


-dS- 


7i(Y,)n • r e"“* e c e 


icar -iaMY j 




-dS (2.21) 


f=0 f=o 

The time derivative of the first term in eqn. (2.21) is calculated at fixed x . It can 
be obtained using the material derivative operator of eqn. (2.11). After algebraic 
manipulations, eqn. (2.21) can be rewritten as a single surface integral in the form 


4rcPP 5 (X,t) = e"'“' e 


- ioMXi 


n (Y,)n i 


ia 


^ 2 + n 2 +? 2 (^+r 1 2 + q 2 ) 3 ' 


exp^ia V^ 2 + Tl 2 + ? 2 )dS ( 2.22) 


f=0 

Equation (2.22) is the integral representation of the scattered pressure written in terms of 
body-fixed coordinates. At this stage, eqn. (2.22) applies for an arbitrary duct cross- 
section. 


In the next phase of the analysis, the incident field is utilized in conjunction with 
eqn. (2.22) to form the boundary integral equation. To accomplish this goal, a boundary 
condition is introduced for a duct with rigid walls. The no penetration condition is 
applied, namely Ut • n = 0 on the duct surface f=0. Substituting the no penetration 
boundary condition into the linearized momentum equation yields 


d(u, : n ) _ dfr = Q ( 2.23) 

dt dn 

Then use of eqn. (2.2) implies that the boundary condition on the scattered pressure is 


= (2.24) 

dn dn 


11 



on the duct surface. Therefore, after taking the normal derivative of eqn. (2.22) and 
utilizing eqn. (2.24), it follows that 


dp - 
dn 


-itiMXi d 

3n I 


«(Y.)n- 


$ 2 + n 2 +<; 2 


(ft’+nVf 


exp^ia-^ + r| 2 + g 2 ^dsi 


( 2 . 25 ) 


in the limit as the observer at X approaches the duct surface f=0. Equation (2.25) is the 
boundary integral equation whose solution for the unknown jump amplitude n is the 
object of the remainder of this work. 

The integral equation (2.25) is now specialized for the rectangular duct by noting 
that, due to the symmetry of the walls and the positioning of the acoustic source, the 
jumps in scattered pressure are identical for Xi = ±a, i.e., the top and bottom surfaces, (see 
Figure 1) and similarly for x 2 = ±b , i.e., the two side surfaces. Thus, the jump amplitude 
7i is defined separately for each independent surface as 


n 


(x,.t) = 


it,(Xj,X 2 ) 

7t 2 (X 3 ,X,) 


on 

on 


X, = ±a 


X 2 = ±bj 


( 2 . 26 ) 


Further, the amplitudes n } and n 2 are represented as expansions in terms of shape 
functions to model the expected oscillatory behavior present on the surface of the duct. 
Thus 7T | and n 2 are written as 


/ x J -' / \ / x 

*i(Y 3, Y 2 ) = X Xa j , k (P j (Y 3 )M/ k (Y 2 ) 

j=0 k -0 ( L.L / ) 

J - ! K2 + K1-1 . . . . 

Jti(Y 3 ,Y,)= I I a j . k «P j (Y 3 )v k (Y 1 ) 

j=0 k = Kl 

where K„ K 2 and J are the number of functions included to describe n on the surfaces 
Xi = ±a and x 2 = ±b and along the length of the duct, respectively. The unknown 
coefficients a j k in eqn. (2.27) are to be determined. T.ie specific forms chosen for the 
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shape functions cp_, and y/ v will be explored in detail in Chapter 3 of this thesis. The 
integral equation for the rigid, rectangular duct then appears in terms of contributions 
from each of its four sides in the form 


5P:(X,t) f o 

a n 1 v ' -l tot-*-- 

= e 1 


J-l K I — 1 r t J-l KI + K2- 1 r -l 

X ZaJlMJ +Z I a j . k [lj?k , + Ijjl] 

1=0 k=0 1 1 j=0 k = KI 1 J 


where 


IJV- T ‘P j (Y 3 )!v k (Y 2 )(X,-a)F(^T 1 ,X l -a)dY J dY3 

-L/2 -b 

Iff- t «P J (Y J )!v k (Y,Xx,-b)F(6.X 1 -b,c)dY l dY J 

-L/2 -• 

1$ = - T fp, (Y,) J (y 2 )(x, + a)F(^ t|, X, + a)dY 2 dY 3 

-L/2 -b 

IJV = - T <Pj(Y,) j M» k (Y t)(x 2 + b)F(^, X 2 + b,c)dY, dY, 

-L/2 -• 


and where the function F is defined by 



In solving eqn. (2.28) later it will be required that the integral equation be 
evaluated for observers on the top surface X,=a and on the side surface X 2 =b. The 
analysis here will be presented only for an observer positioned on the top surface, and it 
will be assumed that X 2 is on the interval -b<X 2 <b so that the observer is never precisely 
at the comer of the duct. The same type of analysis can also be utilized when an observer 
is positioned on the side surface of the rectangular duct. It is unnecessary to repeat the 
derivation however, because the result follows by analogy with that for the observer on 
the top surface by simply interchanging X, and X 2 , Y , and Y 2 , C, and t|, and a and b. 
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The normal (or X,) derivatives indicated in eqn. (2.28) are now calculated to put 
the integral equation in its desired final form. The left-hand side of eqn. (2.28) follows 
from eqn. (2.12) evaluated on X,=a as 


dP.(x,t) 

ax, 


, , — lOUVLA 3 

= H 1 )(x 3 )e — p - e-^' 


(2.31) 


where 


p< 1) (x 3 ) = 


P 0 CA Xl 

^p 2 ( Bo r 



+ ia B 0 



-mx 3 

p 


+ 



( 2.32) 


and B 0 is as previously defined. The integrand in the term 1$ of eqn. (2.28) is singular 
when the observer and source points coincide (i.e., X 3 =Y 3 , X 2 =Y 2 , Xi=Y,=a), but the 


other three terms in eqn. (2.28) are non-singular as long as X 2 ^ b . Since the integrals in 
eqn. (2.28) Eire to be evaluated numerically, the singularity in the first integral must be 
removed and treated analytically. Thus, the normal derivative of 1$ will only be written 
symbolically at this time. On the other hand, the normal derivative of the three non- 
singular terms can be calculated directly for X,=a. The manipulations required are purely 
algebraic and, after dividing out the common exponential factor, the boundary integral 
equation becomes 


where 


-4jiP pf !) = K (l) + K (2) + K (3) +K (4) 


r) J-l K.I -1 J-l 

K (,, = lim — I Ia w I<S=I 

Xi-*a C*X| j=0 k=0 j=0 


K.l -1 

£ 3j M 
k=0 


- f 1 

Xi-»a dX\ 


( 2.33) 


(2.34) 
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d J- 1 

KUK2-1 

J-l KDK2-1 

= lim — — Z 

Z a j,k 

I< 2 k > = -I I a, 

Xi-*« 5Xi j=° 

k = KI 

j=0 k = Kl 

d J-l 

Kl-I 

J-l Kl - 1 L/2 

1! 

!|- 

M 

Zaj,klji 

l“-Z Zaj.k / <( 

5X] 

k=0 

j = 0 k =0 l/2 

f) j-i 

K1 + K2-1 

J-l K1+K2-1 1 

ii 

3 

1 1 
M 

Z a Jt k 

II 

M 

M 

u> 

X|-+a (?X| j*0 

k = Kl 

j=0 k = Kl 1 


(2.35) 


and where the function G is defined by 


G(*L,n.;)= 


3 


3ia 

_ 4_ 

(ia) 2 



5/2 


2 + 

1 

(V+’l 2 -'f) 

3/2 

/ 


(2.36) 


In the following section, the singular term given in eqn. (2.34) will be treated 
analytically. The other three terms involve only non-singular integrals, and these will be 
evaluated numerically in the form in which they are presented in eqn. (2.35). 


2.4 Singularity Analysis 

The singular integral term that appears in eqn. (2.34) contains the integral defined 
in eqn. (2.29). If the symbol h is introduced, where h=X,-a, it can be written as 

lim |-{ h T <P j (y 3 ) J M\( Y 2 )f(^ , t|, h) d Y 2 d Y 3 j < 2 - 37 > 

h-*0 dh h-*() oh [ -U2 -b ) 

Define I k (Y 3 ,h) to be the transverse, inner integral in eqn. (2.37) so that 


l k (Yj,h)= W k (Y 2 )F(^n,h)dY 2 (2.38) 

- b 

This single integral is considered first; it is singular when q and h vanish, which occurs 
when the observer point coincides with the source point on X,=a. By expanding the 
integrand for small values of £,, q and h, the following expression is obtained 


H\F 



Oaf 


A 2 n 

(B) VJ 


- 0 ( 1 ) 


in which B = f + q 2 + h 2 , and 


( 2.39) 
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A 0 = H' k (X 2 ) , A, = v|/' k (x 2 ), A 2 = ^ 

The singular terms shown explicitly in eqn. (2.39) are denoted as (v(/ k F) 0 , and eqn. (2.38) 
is written as 

I t (Yj,h)= j[v k F-(v k F)JdY,+ !(v k F).dY, ( 2 - 40 > 

The first integral in eqn. (2.40) is completely non-singular as £, r\ and h approach 
zero. The second involves only elementary integrations which can be evaluated 
analytically by changing the integration variable to r|=Y 2 -X 2 . This results in 



Here the r\ limits are defined by r| =-b-X 2 and q + =b-X 2 ; therefore, neither limit vanishes if 
-b<X 2 <b and only the last term on the right hand side of eqn. (2.41) is singular as £, and h 
go to zero. 

The expression in eqn. (2.40), including the result of eqn. (2.41), is now 
substituted back into eqn. (2.37). At this stage of the analysis, each of the non-singular 
integrals that appear can be differentiated with respect to h and the limit as h goes to zero 
can be calculated directly. This yields 

= i> i(Y5) !h F -( v, - F ).L dYjdY3 
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( 2 . 42 ) 


dYi 

The first two terms on the right in eqn. (2.42) can be evaluated numerically without 
difficulty, and will not be discussed further. The entire singular behavior of the boundary 
integral equation (2.33) has now been isolated in the third term, and this is the subject of 
the remainder of this section. 

To further analyze the singular term in eqn. (2.42), define the function F 0 as 

( 2 . 43 ) 

Since neither r| + nor q. vanish, the function F 0 is non-singular when ^=h=0. The third 
term on the right side of eqn. (2.42) is rewritten as 



a 

- Aolim— 
h-*o ah 


h f <Pj(Y3) 


($ 2 + h 2 )fe + n 1 + h 1 
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(2.48) 


L \ L r 2 ^(Y 3 ) /„ l y l { 2 4<Pi(Y 3 ) l/ 2 ^ 2 F 2 (^hW(Y 3 ) 

I, = Fo(0, h) J_-^- T7 dY»+Fi(0,h) / - dY 3 + i r^~7 dY 3 


-L/J^ + h 2 


-L/2 % +h 


^ + h 2 


The third integral in eqn. (2.48) is now non-singular at h=0. 
From eqn. (2.46) it follows that 


Fi(0, h) = lim 

i-*o 


F.fe,h)-F„(0.h) 

4 


dF o (0,h) 


(2.49) 


In addition, making use of the fact that p£,=Y 3 -X 3 , it is seen that the first integral in eqn. 


(2.48) can be rewritten as 


It then follows that 


ip <Pj(y 3 ) dY 3 
-l/2 i, 2 + h 2 


hax 3 -L/2 


J 2< Pj( Y 3 ) tan | dY 3 

L/2 vn/ 


( 2.50) 


lim J-(hl s ) = lim 4-(-PFo(0.h)-|- T (p (Y 3 )tan- , f^]dY 3 + h L f - s --?^’ h ^^ Y ^ dY3 | (2.51) 

h-oSh' h-»o5h 5X 3 -l/2 ' VhV - L / 2 ^ +h 


The limits in eqn. (2.51) can now be evaluated. The first term on the right is 


^ < P j (Y 3 )tan-(QdY 3 + PF o (0,h)^- ? |^^dY 3 
h-^0 C/h 5X3 -L/2 y YiJ 0X3 -L/2 C + h 


h->o I dh -L /2 4 + h 2 5 X 3 -l/2 £ + h 2 


(2.52) 


after differentiation. The first term on the right hand side of eqn. (2.52) can yield a non- 
zero limit only because of contributions of the integral at the singular point ^=0. The 
limit is therefore obtained by reducing the integration range to small local interval -5<^<S 
so that it becomes 


lim 

h-+0 




1 


ah 


6^ + h 2 



lim 

h— >0 


'^ l( x,)„-(i)‘ }-o 


(2.53) 
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because dF 0 /dh vanishes at h=0. The second term on the right of eqn. (2.52) is 
expressible in terms of a Cauchy principal value integral as 


r) l/2 <p.(Yi) 

E dYs 

0 X 3 -L/2 C, 

after using the fact that F o (0,0)=2. 

Now, consider the second term on the right side of eqn. (2.51). It is 


'i-f “4- 

h-»0 dh -L/2 q + h | 

12 L/2 

/ F 2 M<P i (Y,)dY3 + lim h / 


L/2 

I 

-L/2 


S 2 

3F 2 (^,h) 

£ 2 + h 2 

ah 


V + h 2 


«Pj(Y,)dY,} 


( 2.54) 


( 2.55) 


It is easily shown, by the same procedure that led to eqn. (2.53), that the local 
contribution of the second term in eqn. (2.55) also vanishes. 

Finally, because 


F„(i;.0)-Fo(0,0)-i;F l (0,0) Fo(U>)-F < ,(0,0) ( 2 , 56 ) 

F 2 (£.°J- ~2 ^2 

it follows from eqn. (2.55) and (2.56) that eqn. (2.44) is 


-A 0 lim— (hl s ) = -A 0 
h->o on 


8 H 2 ^(Y,) L ' 2 

2 P^T-^^ dY 3+ I 

5X 3 -L/2 S -L/2 


<Pj(Y 3 ) 

n 

n+ 

n 

n+ 

: * 2 

Vn 2 + £ 2 

'Ini 

n- 

n- 


dY 3 


( 2.57) 


where r| + and r\_ are as previously defined. The shape functions cpj will be defined in the 
following chapter of this thesis and it will be shown there that an analytical expression for 
the Cauchy principal value integral can be obtained. It should be emphasized here that 
the second integral in eqn. (2.57) is non-singular and can be evaluated numerically 
without any difficulty. 
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2.5 Summary 


At this point it is appropriate to summarize the lengthy analysis just presented. 
Equation (2.33) is the singular integral equation to be solved for the unknowns a^. To 
solve this equation numerically requires the evaluation of the integrands that appear as 
coefficients of the unknowns in eqns. (2.34) and (2.35). However, the singular integral of 
eqn. (2.34) cannot be evaluated directly. The preceding section has outlined the method 
by which the singularity has been removed from the integral. Its complete analytical 
evaluation will be carried out in the following chapter. The numerical task remaining 
involves only the straightforward evaluation of the double integrals in eqn. (2.35), the 
double integral in the first term on the right in eqn. (2.42), the single integral in the 
second term on the right in eqn. (2.42) and the single integral in the second term on the 
right in eqn. (2.57). Again, all of these integrals are non-singular and can be evaluated 
numerically without difficulty. As mentioned previously, all of these expressions can be 
converted to apply to an observer on the side surface by the appropriate interchange of 
variables. 

The next chapter outlines the numerical procedure followed to obtain the 
unknowns aj k which determine the unknown scattered pressure jump across the duct 
walls. Once known they are utilized to calculate the scattered pressure at an arbitrary X 
in the field from 


p,(x,t) 



( 2 . 58 ) 
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where the 1$ are the integrals defined in eqn. (2.29); all of these integrals are non- 
singular so long as X is not on the duct surface. The total acoustic pressure radiated 
from the duct is obtained by addition of the known incident pressure to the calculated 
scattered pressure. 
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3. Numerical Implementation 

The previous chapters detailed the development of the set of equations which, 
when solved, predict the acoustic field generated by a monopole source placed at the 
center of a rectangular duct. The equations have been analytically treated to eliminate 
problems in implementing numerical schemes. The governing equations for this problem 
have been coded into a FORTRAN code and run on a 500MHz DEC-Alpha workstation. 

Before addressing the issues in this chapter it is pointed out that the axial 
integration in eqn. (2.33) is similar to that done by Kosanchick [4] in applying a boundary 
integral technique to predict the acoustic field generated by propellers in lined as well as 
rigid circular ducts. A shape function expansion was also posed there for the axial 
behavior of the scattered pressure jump along the surface of the duct. However, the form 
of the expansion proved to require adaptive techniques at the leading and trailing edges of 
the duct due to sensitivity to the choice of parameters. The code developed for the rigid 
circular cylinder was later modified to replace the propeller noise source model with a 
point monopole. Results describing radiation from the monopole acoustic source were 
discussed by Myers [6]. Prior to the current consideration of the rectangular duct, a new 
set of shape functions were introduced to describe scattered pressure jump along the 
surface of the circular duct. Through extensive numerical testing, it was found that this 
new set of shape functions eliminated the need for additional adaptive techniques at the 
leading and trailing edges of the duct. The current analysis for the rectangular duct 
makes use of this new set of axial shape functions and they are discussed in the following 
section of this chapter. 
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3.1 Pressure Jump Representation 


The scattered pressure jumps are represented as a series of products of two 
sequences of carefully chosen shape functions as presented in eqn.(2.27). The shape 
function expansions are chosen to model the oscillatory behavior expected on the duct 
surface as well as the appropriate behavior of the pressure jump at the leading and trailing 
edges of the duct. With these criteria in mind, the shape function expansion for the axial 
direction is discussed first. 

It is known that integral equations like the one under consideration here do not 
have unique solutions until further conditions associated with the edge behavior of the 
solution are specified. As in thin airfoil theory, the pressure jump that satisfies eqn. 
(2.33) can be expected to have an inverse square-root singularity at the leading edge of 
the duct when M * 0 , and a square-root zero is anticipated at the trailing edge of the duct 
[11]. This latter condition is the well known Kutta condition. These edge conditions can 
best be imposed by introducing the variable 

Y 3 = — ~cosk (3.1) 

where 0 < k < n , and by defining 

X 3 = --^cosk„ ( 3 - 2 ) 

Now, the axial variation of the scattered pressure jump is written in terms of the set of 

functions commonly utilized in thin airfoil theory [1 1]: 

l + COSK . 

cP-(k) = • sinK J (3.3) 

sin jK, j = l,...,J-l 
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The leading edge singularity at k=0 is contained explicitly in the first term of the 
sequence. It is easily shown using eqn. (3.1) that it grows as the inverse square root of 
distance from the edge. The remaining terms vanish at the leading edge of the duct and 
all vanish as the square root of distance from the trailing edge in accordance with the 
Kutta condition. As mentioned previously, use of eqn. (3.3) with a uniform discretization 
in k rather than Y 3 , has been shown to eliminate the need for additional adaptive 
techniques in handling the square root singularity at the leading edge of the duct. Most 
importantly, the Cauchy principal value integral of eqn. (2.54) can be obtained 
analytically for the shape function expansion posed in eqn. (3.3). This analytical result 
will be discussed in the next section. 

The variation in the scattered pressure jump in the lateral direction on the wall 
surfaces is expected to be similar to that seen in the incident pressure field. This variation 
can be replicated by a sequence of functions that are sinusoidal in form. Therefore, the 
shape functions modeling the lateral oscillations of the jump on the surfaces of the duct 
are taken in the form 


and 


v k (v 2 ) 



k =0 
k =1 


( 3 . 4 ) 


II 

> 

l 

\ x f Yi + a V 


sin 

^(2k-2K,-l)«[ 2a Jj 


k =Ki 

k = K.i + 1 Ki + K.2 - 1 


( 3 . 5 ) 


where K, and K 2 are the finite number of functions to be used in the sequences. Both sets 
of functions are symmetric about the midpoints of their respective lateral surfaces. 
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It is noted that only a finite number of functions are used. Choosing the correct 
number of functions from each shape function expansion set is not a straightforward task. 
The number of functions that is sufficient for an accurate solution varies depending on the 
choice of problem parameters such as Mach number, frequency and length of duct in 
relation to the cross-sectional geometry. The process of determining the proper number 
of functions is discussed later in conjunction with a numerical discussion of a sample 
problem. 

3.2 Principal Value Integral 


The Cauchy principal value that appears in eqn. (2.57) can be evaluated 
analytically for the shape functions defined in eqn. (2.33). By utilizing the relation 
P^=Y 3 -X 3 and the transformation of eqn. (3.1), the integral becomes 


F) L/2 

|( | pv) = 2p 2 — — X 


dXi -l/2 (Yj-Xj) 




0X30 COS K “ COSKo 


-sin KdK 


(3.6) 


The shape functions are substituted into eqn. (3.6) and the known integral 


" cosnK , nsinntCo 
J( dK = 


0 COSK -COSKo 


SUlKo 


derived, for example, by Karamcheti [11], is used. For j-0 this yields 


Io PV) = -2p 2 — — rf / 1 + C0SK = -2[f— ^-[0 + 7l] = 0 


9X30 (cos K -COSKo) 


9Xi 


When j>0, eqn. (3.6) is 


](P v ) = -2p 2 


9 ^ sinjKsinK 
9X3 0 (cos tc - COS K 



Utilizing the trigonometric identity 


(3.7) 


(3.8) 


(3.9) 
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( 3 . 10 ) 


sin jk sin k = 


cos( j - 1 )k - cos( j + 1 )k 


and the relation established in eqn. (3.7), eqn. (3.9) is 


ir>=2 P ^( 7I co S jQ= - 47tp2jsinjK ° 

J H 5X 3 V ’ L sin K „ 


( 3 . 11 ) 


where eqn. (3.2) has been utilized to obtain the derivative with respect to X 3 . The results 
expressed in eqns. (3.8) and (3.1 1) constitute, in analytical form, the entire effect of the 
singularities in the integral equation (2.33). 


3.3 Method of Collocation 


Given a choice of the number of shape functions J, K,, and K 2 , there are a total of 
J(K,+K 2 ) unknown coefficients a j k to be determined to complete the solution of eqn. 
(2.33). Although there are numerous techniques available, the method of collocation [12] 
is the technique utilized here to obtain the solution. This method involves the selection of 
a grid of observer points on the surface at which the integral equation eqn. (2.33) is 
satisfied exactly, thus producing a set of inhomogeneous linear algebraic equations for 
the unknown coefficients. In the current work, however, the method is supplemented by 
imposing the obvious physical condition that the pressure jumps on each wall should 
equal one another at the comers where the walls meet. Thus, for example, it will be 
required that 7ti(X3,b) = 7t 2 (X 3 ,a) . If an unlimited number of functions were used in 
each series expansion, this condition would presumably not be necessary. However, the 
constraint imposed by the comer condition has been found to lead to an accurate 
representation of the jump in scattered pressure using what appears to be as few shape 
functions as possible. 
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The grid of observer points used here consists of J locations evenly distributed 
along the X 3 direction of the surface over the interval -L/2< X 3 <L/2. Because of the 
symmetry present in the current problem, only one quarter of the perimeter of the duct 
must be considered. Thus, K, evenly distributed observer locations are chosen over the 
interval 0 < X 2 < b on the top surface, and K 2 -l locations are distributed over 0 < Xi< a 
on the side. The integral equation is satisfied at these J(K,+K 2 -1) points. The remaining 
J equations necessary to complete the algebraic system are obtained by imposing the 
comer condition at the J axial stations of the grid. 

The first and last axial collocation locations are at a distance D* in from the 
leading and trailing edges of the duct. The first grid point for each side is placed at the 
center of each side (i.e., X 2 =0, X,=0) and the last is chosen at a distance of 10% of the 
half length of the side in from the comer. 

3.4 Numerical Integration 

With the set of collocation points defined, the numerical evaluations of the 
integrals of eqn. (2.34) is discussed. Four point Gauss-Legendre Quadrature is utilized to 
complete all of the necessary integrations [12,13]. The integrands in eqn. (2.34) are all 
non-singular and are well behaved and can be numerically integrated accurately given 
that the proper discretization is utilized. The discretizations required along Y 2 and Y, are 
independently determined based on the oscillatory behavior of the integrands. For the 
surfaces Xi = ±a , the oscillations of the integrand are governed by the shape functions 
v}/ k (Y 2 ) and the exponential function in terms of £, and r| as 
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sin 


(2k - 1 )rc 


Y; + b 
2b 




( 3 . 12 ) 


The Y 2 oscillations are most rapid when £,=0 and Q= 0. If the sine term in eqn. (3.12) is 
rewritten in exponential form, the fastest oscillatory behavior with respect to Y 2 is 


; (i(;>k-i)*(^) + Mv.-xj) = e (^k( 


(2k-l)n 


2b 

Therefore, the wavelength of the fastest oscillation is 


~aX : 


( 3 . 13 ) 




2n 


(2k -1) 
a + — — — 7t 


2b 


( 3 . 14 ) 


The total number of panels along the X t = ±a surface is calculated using 



( 3 . 15 ) 


where P 3 is the number panels chosen per wavelength. This parameter allows for the 
optimization of the level of discretization utilized. Two panels per wavelength were 
chosen. The same procedure determines the most rapid variation in the oscillatory 
behavior of the integrand in eqn. (2.34) on the surfaces X 2 = ± b . 

The axial discretization for the Y, integration of all of the functions expressed in 
eqn. (2.33) is examined now. The relationship between the variables Y 3 and k was 
established in eqn. (3.1). Utilizing an even distribution of integration points in k leads to 
a clustering of integration points in Y 3 around the leading and trailing edges of the duct 
eliminating the need for additional adaptive techniques. The length of the duct is 
discretized evenly in k on the interval [0 ,ti]. The variation in the axial direction is 
described by 
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(3.16) 


sin jKe ia ^ 2+r,+ ^ 2 

The oscillations in k are most rapid when t|=0 and ^=0. If the sine term in eqn. (3.16) is 
written in exponential form and § is written in terms of k, then the most rapid variation in 
the axial direction is approximated by 


Now the wavelength 


(L ) 

e i£ * e ijK »e i ll a+J J K 


of the fastest oscillation is determined by 

2n 




L 

a — + j 
2 


(3.17) 


(3.18) 


The total number of panels required for the integration in the axial direction is calculated 
using 


N k = P>.— ( 3-19) 

where again, P x is the number of panels per wavelength. Four panels per wavelength are 
used for all of the axial integrations. The integration with respect to Y 3 is rediscretized for 
every j and, likewise, the integration with respect to either Y, or Y 2 is rediscretized for 
every k. The goal is to utilize as coarse a discretization as possible while maintaining a 
uniform level of accuracy. The testing of the accuracy of the level of integration utilized 
will be explored in the following section. 

At this point, the numerical techniques necessary to obtain the solution to the 
integral equation in eqn. (2.33) have been described. In the next section the solution of a 
sample problem is investigated. 
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3.5 Sample Problem Solution 


The previous sections outlined the methodology necessary to determine the 
acoustic field generated by a monopole source placed in a rectangular duct. The 
methodology is applied here to a sample problem to test the integrity of the FORTRAN 
code. The following table indicates the parameters chosen for the sample problem: 


Parameter 

Value 

Parameter 

Value 

M 

0.0 

a 

0.25 m 

L 

0.5 m 

b 

0.5 m 

Frequency 

750 Hz 

c 

340. 1 7 m/sec 


Here, the frequency and duct dimensions are relatively small in order to minimize the 
numerical effort required for the testing, and the duct and source are stationary in the 
fixed medium. In this case, the incident field is symmetric in the axial direction about 
X 3 =0 and therefore the jump in scattered pressure is also expected to be symmetric about 
X 3 =0. Also, when M=0.0 there is no singularity in scattered pressure at the leading edge 
of the duct. Thus, for this test case, the coefficients of the singular first terms and the 
terms that are anti-symmetric about X 3 =0 (j=2,4,6...) in the axial shape function 
expansions (2.27) must be driven to zero in the solution of the algebraic system that 
results from collocation. 

The incident monopole pressure field is illustrated in Figures 2-4. Figure 2 shows 
the axial variation in the complex incident pressure amplitude in both real and imaginary 
parts on the duct surface X,=a along its centerline X 2 =0. The complex incident pressure 
amplitude for the lateral centerline X 3 =0 is shown in Figure 3 for the surface X,=a and is 
shown in Figure 4 for the surface X 2 =b . 
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The method of collocation is now used to solve for the coefficients % k in eqn. 
(2.33). The grid of collocation locations is as discussed in section 3.3, and the axial 
offset D* is taken to be 0.1L. Assuming that the collocation grid is adequate and that the 
discretization discussed in section 3.4 yields accurate integrals, the number of functions 
from the expansions (2.27) required to ensure a converged solution must be determined. 
The oscillatory nature of the incident data function Pi (l) (X 3 ) in eqn. (2.33) dictates the 
number of functions required. However, this data function oscillates in a manner similar 
to that of the incident pressure so that the number of functions required can be inferred 
from the oscillations in the real and imaginary parts of the incident pressure illustrated in 
Figures 2-4. It is seen in Figure 2 that about 1 'A oscillations exist in the incident pressure 
along the length of the duct for the sample problem. Experience indicates that it is 
necessary to include all terms in the expansions (2.27) up to terms that oscillate at least 
twice as fast as the incident data [4]. The axial shape functions (3.3) for j > 1 each have j 
half-oscillations over the length of the duct. Therefore, a minimum of 6 of these sine 
functions are expected to be required to reproduce the 1 l A wavelengths (J=7). A total of 
8 functions are actually chosen. Similar considerations determine K, and K 2 from the 
incident pressure amplitude plotted along the lateral centerlines X 3 =0 on the surfaces 
X,=a and X 2 =b. These lead to the choices K,=8 and K 2 =6. Upon solving eqn. (2.33), the 
coefficients a, k are obtained from which the scattered pressure jump is calculated utilizing 
eqn. (2.27). The resulting pressure jump amplitude along X 2 =0 on the surface X|— a is 

shown in Figure 5. It is seen that, as expected, the jump amplitude is symmetric about 
X 3 =0 and there is no singularity present at the leading edge of the duct. Thus the scheme 
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has driven the singular terms and the anti-symmetric terms in the expansions (2.27) to 
zero as required. The jump amplitude n ] on the surface X,=a laterally along X 3 =0 is 
shown in Figure 6 and the jump amplitude n 2 on the surface X 2 =b along X 3 =0 is shown in 
Figure 7. In Figures 5,6 and 7 the jump amplitudes are all smooth and, as expected, 
display oscillatory behavior similar to that of the corresponding plot of the incident field. 

To determine if the solution obtained is converged, the number of functions used 
in the expansions (2.27) is varied. First, K, and K 2 are held constant and J is varied as 
shown in Figure 8. As J is increased there is little change in the scattered pressure jump. 
A conclusion can be made here that a sufficient number of axial shape functions were 
chosen. The same procedure is followed for the lateral functional expansions. Figures 9 
and 10 show only minor changes when holding J=8 and varying the values of K, and K 2 
for the X 3 =0 lateral plane on the surface X,=a and on the surface X 2 =b, respectively. 
Therefore, it can be concluded that using J=8, K,=8 and K 2 =6 is sufficient to accurately 
represent the solution for the jump amplitude in this case. Further tests involving 
increasing all three parameters simultaneously were carried out, and all supported the 
same conclusion. 

It is also necessary to verify the validity of the discretization utilized. This can be 
achieved by either increasing the number of Gauss points or by increasing the number of 
panels per wavelength, P x . Here the number of Gauss points is doubled from 4 points to 
8 points. Figures 1 1 and 12 show that no apparent differences exist between results with 
the 4 point and 8 point Gauss-Legendre schemes along the lines X 2 =0 and X 3 =0 on the 
surface X,=a. Many other such plots were examined, and all indicated that the level of 
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discretization utilized was sufficient to obtain accurate results. The remainder of the 
results presented in this thesis are obtained utilizing the 4 point Gauss-Legendre scheme 
with P x =2 for discretizing X, and X 2 and P x =4 when discretizing in X 3 . 

Finally, the choice of collocation points is considered. The collocation points are 
evenly distributed along the length of the duct and along half of each lateral side. The 
axial locations are distributed evenly between points at a fixed distance D* in from each 
end of the duct. A different set of collocation points results simply by changing the value 
of D*. Figure 13 shows that the solution with D*=0.005L exhibits no significant 
variations from the solution with D*=0.01L along the axial centerline of the surface X,=a. 
Other points with varying D* values showed little variation in the solution as well. 
Similarly, the first lateral collocation points are on the center of each side and the last is 
placed at a distance 10% of the half-length of the respective side. The set of collocation 
points is modified by changing the position of the last point. The last collocation point 
was placed at a distance of either 9% or 1 1 % of the half-length of the lateral surfaces. 
Figure 14 shows only slight changes in the solution when the set of collocation points is 
altered. The distance corresponding to 10% of the half-length is deemed reasonable since 
variations as small as these seen in Figure 14 have an insignificant effect on the radiated 
field which is the topic of the next chapter. 

It should also be mentioned here that if the dimensions a and b were taken to be 
equal then, upon solving the system of equations, the jump amplitudes 7t, and n 2 should 
be identical. When the parameters selected for the sample problem were utilized in 
conjunction with the condition that the dimension a was the same as dimension b, the 
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jump amplitudes were found to be the same as expected. The same result was obtained 
when the observer points were taken on one surface and thereby reducing the number of 
unknowns. 

3.6 Moving Duct Problem 

While the problem discussed in the previous section was a useful test case, it is 
also of interest to examine the solution for the jump amplitude when the duct is in motion 
so that the jump is singular at the duct leading edge. Therefore, a case for which M=0.1 
is presented here. The same numerical tests were also carried out for this case and led to 
similar conclusions in regard to the accuracy of the level of discretization as well as the 
choice of collocation locations. They will not be discussed further. For the moving case, 
the same parameters as utilized in the previous problem were retained except that the duct 
length was taken to be 2m. This is four times the length used in the sample problem and 
therefore 4 times as many sine functions in the axial shape function expansion would be 
expected to be required to obtain an accurate solution. However, for this flow case 25 
functions from the expansion (2.27) were not sufficient. In fact, it was found that a total 
of 52 axial functions were required before the solution was completely converged. It 
appears that the need for so many functions arises because, in contrast to the circular duct 
treated in ref. [6], the singularity at the leading edge varies in strength around the 
perimeter of the duct. This differs significantly from the circular case where the 
singularity was determined by a single point at the leading edge of the duct. This is one 
reason why the rectangular duct problem is much more computationally intensive than 
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the circular case treated in ref. [6]. A total of K,=8 and K 2 =6 functions from the lateral 
surface expansions were found to be sufficient and no apparent change in the solution 
was observed when these values were varied. Figure 15 shows the scattered pressure 
jump on the surface X,=a along X 2 =0 for a rectangular duct in the presence of flow. 
Notice the square-root singularity at the leading edge of the duct and the square-root zero 
at the trailing edge of the duct. 

The radiated field for the above sample problems, as was well as for a higher 
frequency case are presented in the next chapter. 
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4. Radiated Field 


The previous chapters outlined the analysis necessary to determine the acoustic 
field generated by a monopole source placed at the center of a rectangular duct and some 
of the numerical checks made on the analysis for inconsistencies. Once the jump 
amplitudes rc, and n 2 have been determined, the scattered pressure at any observer 
location is obtained by substituting these jump amplitudes into eqn. (2.22). The total 
pressure at any point in the field is then obtained through the simple addition of the 
incident and scattered pressure at an observer point in the field. The actual or physical 
pressure for the incident, scattered or total pressure is the real part of its complex 
counterpart. The corresponding root-mean squared (RMS) pressure is calculated in the 
usual manner such that 


fpp 7 

Pr ms = ]lY' < 41 > 

where p* is the complex conjugate of p. The sound pressure level (SPL) is determined 
by 


SPL = 20 log 


Pnns 
' Pref ' 


( 4 . 2 ) 


where p ref =2 x 1 O' 5 Pa. 

The incident, scattered and total pressures are illustrated in the following results 
one cross-sectional plane at a time on the intersection of a sphere centered at the origin of 
the body-fixed coordinate system and the plane of interest. These intersections are circles 
along which are give:: polar plots or directivities of the radiated field about the plane of 
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interest in terms of sound pressure levels in decibels (dB) from eqn. (4.2). A radius R of 
5m was utilized for all of the directivities presented in this thesis. Figures 16(a)-(c) 
depict the center planes of interest. The directivity angle is called 0 in each case, and it is 
measured from the positive X„ axis (O°<0<18O°). This angle is utilized to obtain 
directivity plots for planes that intersect the source so that the radius of the arc for these 
polar plots is 5m. For planes that do not intersect the source, such as that shown in 
Figure 16(d), the radius of the polar arc is less than 5m. Since the incident field is 
spherically symmetric when M=0, repeated use of the same spherical radius leads to an 
incident field that is the same for each directivity plot for this case. 

Chapter 3 included a detailed numerical validation process as well as the 
presentation of the jump amplitudes for the sample problem. The radiated field plots for 
the sample problem in the center planes X 2 =0, X 3 =0, and X,=0 are presented as solid lines 
in Figures 17 (plane in Figure 16(a)), 18 (plane in Figure 16(b)) and 19 (plane in Figure 
16(c)), respectively. The symmetric nature of these field plots is expected due to the 
symmetry in the pressure jumps seen earlier for this case. For a duct length of 0.5m, the 
solid lines in Figures 17 and 18 illustrate significant differences between the incident and 
scattered pressures at 0=90°. These give rise to a substantial enhancement in the total 
field about 0=90° in the plane normal to the wide side (normal to surface X,=a) of the 
duct. There is less apparent difference between the incident and scattered pressure in 
Figure 19 resulting in only a slight enhancement about 0=90° in the plane normal to the 
narrow side (normal to surface X 2 = b) of the duct. Figures 17 and 19 also indicate that the 
monopole field is nearly unaffected by the duct in the direction of the central X 3 axis. 
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The enhancement seen in Figure 17 is unexpected since the results for the circular 
ducted fan problem presented in refs. [2-4] usually depicted shielding rather than 
enhancement about 0=90°. However, the propeller in ffee-space radiates primarily in the 
lateral direction so that its radiation is significantly altered by the presence of the duct. 
The monopole radiates in all directions including axially to the ends of the duct. This 
difference in source type could explain the enhancement shown in Figures 17-19. In Ref. 
[6], the radiated field from a monopole in a circular duct was discussed. The code 
utilized to obtain those results was modified as discussed in Chapter 3 and an identical 
source to that utilized in obtaining the results for Figures 17-19 was specified. Figure 20 
illustrates the directivity pattern for the scattered and total pressure fields for the circular 
duct at M=0. It depicts a small enhancement in the total pressure field about 0=90° for a 
duct length of 0.5m, and indicates that lateral shielding is obtained only if the duct length 
is increased to lm or more. Thus, monopole radiation from a circular duct can also give 
an enhancement like that shown in Figure 17. Furthermore, it appears that the potential 
for shielding by the rectangular duct does exist but that the duct length needs to be 
increased. 

The broken lines in Figures 17, 18 and 19 illustrate the effect on the scattered and 
total pressure directivities of increasing the rectangular duct length to lm and 2m. As the 
duct length is increased, there is less apparent difference between the incident and 
scattered fields which indicates that the potential exists for phase cancellation in the total 
field. It is seen that some shielding is now obtained over the range 3O°<0<15O°. 
Shielding is observed in the total pressure mainly in the center plane X,=0 (normal to 
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narrow side) in Figure 19 and little shielding is observed in the center plane X 2 =0 (normal 
to wide side) in Figure 17. It also appears that shielding exists only for a duct length of 
2m or more. This is not surprising since the center plane normal to the narrow side 
exhibits less enhancement in the 0.5m case. It might be expected intuitively that the wide 
side of the duct should provide greater shielding from the source. However, this example 
shows this not to be true in general. In fact, experience indicates that few general 
statements can be made about the radiation patterns from ducts at low to moderate 
frequencies. Figures 17 and 19 also show that the duct has little apparent effect in the 
more nearly axial directions in the regions 15O°<0<18O° and 0°<9<30°. 

Figures 21 and 22 depict the radiated field plots for various axial and lateral cross- 
sectional planes for the moving duct problem as discussed in Chapter 3. The axial center 
planes of Figure 21 (planes in Figures 16(a) and (c)) and lateral center plane X 3 =0 of 
Figure 22 (plane in Figure 16(b)) intersect the source. Some shielding is obtained in 
these planes in the total pressure about 0=90°. The moving duct does not appear to have 
as pronounced a channeling effect on the total radiated field as was the case for the 
circular duct propfan model in refs. [2-4]. This is due to the difference in source type. 
However, Figure 21 depicts two distinct bulges that display the impact of the singularity 
of the pressure jump on the leading edge and satisfaction of the Kutta condition at the 
trailing edge of the duct when M * 0. These are no doubt counterparts of the channeling 
effect seen for the propfan, but they only exist here in the plane normal to the narrow side 
of the duct. The figure also shows that the level of shielding provided by the wide and 
narrow sides of the duct does not differ as much in the presence of flow as compared to 
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the differences witnessed in the M=0 case in Figures 17 and 19. Because of the 
singularity at the leading edge of the moving duct, it is of interest to examine the 
directivity in the plane that contains the leading edge (X 3 =-L/2). This plane is shown in 
Figure 16(d). The directivity in the plane containing the leading edge along with the 
corresponding plane containing the trailing edge (X 3 =L/2) is shown in Figure 22. The 
incident field shown in Figure 22 is that for the X 3 =0 plane since the incident fields for 
the leading edge, trailing edge and lateral center planes differ only slightly at M=0.1 . For 
a duct moving in the negative X 3 -direction, there appears to be greater lateral shielding at 
the trailing edge than at the leading edge plane, but in general there are not many 
differences in the radiation patterns in these planes. 

Finally, to include a case more nearly representative of practical inlet 
experiments, the frequency is increased to 1922 Hz. At this frequency, the incident field 
is highly oscillatory, and a significant number of shape functions is required from the 
expansions (2.27) to model this field accurately (at least K,=20, K 2 =10, J=12 as a 
minimum). Limitations on time and computational resources lead to the decision to 
reduce the dimensions of the duct and thereby reduce the number of functions required to 
accurately model the field. The dimensions a and b were reduced by a factor of 4 
(a=0.0625m, b=0. 125m, L=0.5m) and only the stationary case M=0 was considered. This 
allowed for accurate results to be obtained with K,=8, K 2 =6 and J=12. Figure 23 (planes 
in Figures 16(a) and (c)) predicts significant shielding in the axial center plane normal to 
the wide side of the duct as would be expected because the effect of the duct at high 
frequencies is generally to beam sound in the axial directions. It is also of interest that, 
despite the relatively high frequency, virtually no shielding is provided by the narrow 


40 



side of the duct, so that the beaming effect occurs here only in planes normal to the wide 
side. Note that the beaming occurs even though the duct is of very short length. In the 
lateral planes, there is significant shielding in the region from 12O°<0<6O° as shown in 
Figure 24 (planes in Figures 16(b) and (d)). Thus the beaming seen in the plane X 2 =0 in 
Figure 23 extends over a rather large angle range in the directions above the wide side of 
the duct. 
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5. Concluding Remarks 


This thesis has detailed the analysis necessary to model the radiated field 
generated by a monopole source placed at the center of a moving rectangular duct 
utilizing a boundary integral technique. The total pressure in the acoustic field was 
written in terms of incident and scattered pressure components such that the scattered 
pressure represents modifications to the incident pressure due to the presence of the duct. 
The scattered pressure is discontinuous across the duct walls. It satisfies a generalized 
wave equation with a source term involving the unknown pressure jump across the duct 
walls. Through the use of an integral representation for the scattered pressure, a singular 
boundary integral equation governing the unknown jump in scattered pressure is derived. 
The jump is represented by a double series of carefully chosen shape functions, and a 
solution to the integral equation is obtained using the method of collocation. With the 
jumps determined, their substitution back into the integral representation of the pressure 
produces the radiated pressure at any observer location. 

The governing equations for the analysis presented in Chapter 2 were 
programmed in a FORTRAN code which was run on a 500MHz DEC-Alpha workstation. 
The singular analysis presented in Chapter 2 isolated the effect of the singularity in the 
set of boundary integral equation thus eliminating problems with the implementation of 
numerical schemes. Numerical verification of the FORTRAN code was presented in 
Chapter 3. The chosen sets of shape functions were discussed and it was shown that 
these functions allowed for a complete analytical evaluation of the singular part of the 
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boundary integral equation. All the necessary numerical integrations were completed on 
non-singular integrals. 

Radiated field results for the some problems of interest were presented in Chapter 
4. For a stationary duct at a frequency of 750Hz, lateral shielding was obtained for a duct 
of sufficient length and was predominantly in the plane normal to the narrow side of the 
duct. In fact, very little shielding was obtained in the plane normal to the wide side of the 
duct, a result which is perhaps not expected intuitively. However, for a higher frequency 
case of 1922Hz, significant shielding by the stationary duct was obtained in the plane 
normal to the wide side of the duct. Interestingly enough, no shielding was provided by 
the narrow side of the duct at this relatively high frequency. 

For a moving duct, the radiated field for the monopole source can be compared to 
that of the propeller source. The channeling effect witnessed with the moving ducted 
propfan was not immediately apparent in the moving ducted monopole case except for 
the appearance of two distinct bulges in the radiated field in the plane normal to the 
narrow side of the duct. At a frequency of 750Hz, comparable shielding was provided by 
all sides of the moving duct. 

Because the analysis required to determine radiation from a rectangular duct is 
fully two-dimensional, it requires significantly greater computational resources than for 
the circular duct. Time and computer resources have allowed for only a few 
representative examples to be studied in this thesis. They almost certainly do not fully 
depict all of the interesting phenomena associated with radiation from a rectangular duct. 
A true understanding of the radiated field would require a detailed parametric study 
involving variations in duct dimensions, frequency and Mach number and it would be 
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beneficial to complete such a study. Furthermore, it would be a logical further step to 
include an acoustic liner in the rectangular duct. Unfortunately, computational 
requirements for this problem would greatly exceed those for the rigid walled duct. 
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Axial Location X3/L 


Figure 2: Incident Pressure Amplitude in plane 
X 2 =0 on Duct Surface X,=a; a=0.25m, 
b=0.5m, L=0.5m, f=750Hz, M=0.0 
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Lateral Location Xj/b 

Figure 3: Incident Pressure Amplitude in plane 
X 3 =0 on Duct Surface X t =a; a=0.25m 
b=0.5m, L=0.5m, f=750Hz, M=0.0 
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Figure 4: Incident Pressure Amplitude in plane 
X 3 =0 on Duct Surface X 2 =b; a=0.25m 
b=0.5m, L=0.5m, f=750Hz, m=0.0 
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Figure 5: Scattered Pressure Jump Amplitude in 
plane X 2 =0 on Duct Surface X,=a; a=0.25m 
b=0.5m, L=0.5m, f=750Hz, M=0.0 
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Lateral Location XJb 

Figure 6: Scattered Pressure Jump Amplitude in 
plane X 3 =0 on Duct Surface X,=a; a=0.25m 
b=0.5m, L=0.5m, f=750Hz, M=0.0 
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Figure 7: Scattered Pressure Jump Amplitude in 
plane X 3 =0 on Duct Surface X 2 =b; a=0.25m, 
b=0.5m, L=0.5m, f=750Hz, M=0.0 


53 






Axial Location Xj/L 


Figure 8: Scattered Pressure Jump Amplitude in 
plane X 2 =0 on Duct Surface X,=a; a=0.25m 
b=0.5m, L=0.5m, f=750Hz, M=0.0, K 1= 8, K 2 =6 
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Lateral Location Xg/b 


Figure 9: Scattered Pressure Jump Amplitude in 
plane X 3 =0 on Duct Surface X,=a; a=0.25m 
b=0.5m, L=0.5m, f=750Hz, M=0.0, J=8 
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Lateral Location X/a 


Figure 10: Scattered Pressure Jump Amplitude in 
plane X 3 =0 on Duct Surface X^b; a=0.25m 
b=0.5m, L=0.5m, f=750Hz, M=0.0, J=8 


56 




Figure 1 1 : Real Part of Scattered Pressure Jump 
Amplitude in plane X 2 =0 on Duct Surface X,=a 
using Different Order Gauss-Legendre Schemes 
a=0.25m, b=0.5m, L=0.5m, f=750Hz, M=0.0 



Figure 12: Real Part of Scattered Pressure Jump 
Amplitude in plane X 3 =0 on Duct Surface X,=a 
using Different Order Gauss-Legendre Schemes 
a=0.25m, b=0.5m, L=0.5m, f=750Hz, M=0.0 
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Figure 13: Real Part of Scattered Pressure Jump 
Amplitude in plane X 2 =0 on Duct Surface X,=a 
with Different Sets of Collocation Points; 
a=0.25m, b=0.5m, L=0.5m, f=750Hz, M=0.0 



Figure 14: Real Part of Scattered Pressure Jump 

Amplitude in plane X 3 =0 on Duct Surface X^a 
with Different Sets of Collocation Points; 
a=0.25m, b=0.5m, L=0.5m, f=750Hz, M=0.0 
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Axial Location X/L 

Figure 15: Scattered Pressure Jump Amplitude in 
plane X 2 =0 on Duct Surface X^a; a=0.25m, 
b=0.5m, L=2.0m, f=750Hz, M=0.1 
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(a) Side View 


(b) End View 


4 L > 



Figures 16 (a)-(c): Spherical Radius for Radiated Field 



< L/2 > 

Figure 16 (d): Intersection of Plane and Sphere 




Figure 17 (a): Scattered Field in plane X 2 =0 at 
Spherical Radius of 5m; a=0.25m, 
b=0.5m, f=750Hz, M=0.0 



Figure 1 7 (b): Total Field in plane X^O at 
Spherical Radius 5m; a=0.25m, 
b=0.5m, f=750Hz, M=0.0 
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Figure 18 (a): Scattered Field in plane X 3 =0 at 
Spherical Radius 5m; a=0.25m, 
b=0.5m, f=750Hz, M=0.0 



Figure 18 (b): Total Field in plane X 3 =0 at 
Spherical Radius 5m; a=0.25m, 
b=0.5m, f=750Hz, M=0.0 
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Figure 19 (a): Scattered Field in plane X^O at 
Spherical Radius 5m; a=0.25m, 
b=0.5m, f=750Hz, M=0.0 



Figure 19 (b): Total Field in plane X^O at 
Spherical Radius 5m; a=0.25m, 
b=0.5m, f=750Hz, M=0.0 
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Figure 20 (a): Scattered Field in plane X 2 =0 at 
Spherical Radius 5m; M=0.0, 
Circular Duct Radius=1.0m 



Figure 20 (b): Total Field in plane X 2 =0 at 
Spherical Radius 5m; M=0.0, 
Circular Duct Radius=1.0m 
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Figure 21 (a): Scattered Field at Spherical Radius 
5m for Moving Duct; a=0.25m, 
b=0.5m, L=2.0m, f=750Hz, M=0.1 



Figure 21 (b): Total Field at Spherical Radius 
5m for Moving Duct; a=0.25m, 
b=0.5m, L=2.0m, f=750Hz, M=0.1 





Figure 22 (a): Scattered Field at Spherical Radius 
5m for Moving Duct; a=0.25m, 
b=0.5m, L=2.0m, f=750Hz, M=0.1 



Figure 22 (b): Total Field at Spherical Radius 
5m for Moving Duct; a=0.25m, 
b=0.5m, L=2.0m, f=750Hz, M=0.1 
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Figure 23 (a): Scattered Field at Spherical Radius 
5m; a=0.0625m, b=0.125m, 

L=0.5m, f=1922Hz, M=0.0 



Figure 23 (b): Total Field at Spherical Radius 
5m; a=0.0625m, b=0.125m, 
L=0.5m, f=1922Hz, M=0.0 
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Figure 24 (a): Scattered Field at Spherical Radius 
5m; a=0.0625m, b=0.125m, 

L=0.5m, f=1 922Hz, M=0.0 



Figure 24 (b): Total Field at Spherical Radius 
5m; a=0.0625m, b=0.125m 
L=0.5m, f=1 922Hz, M=0.0 
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